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Abstract 

A fluid of N smooth, hard spheres is considered as a model for normal (elastic collisions) and 
granular (inelastic collisions) fluids. The potential energy is discontinuous for hard spheres so the 
pairwise forces are singular and the usual forms of Newtonian and Hamiltonian mechanics do not 
apply. Nevertheless, particle trajectories in the N particle phase space are well defined and the 
generators for these trajectories can be identified. The first part of this presentation is a review of 
the generators for the dynamics of observables and probability densities. The new results presented 
in the second part refer to applications of these generators to the Liouville dynamics for granular 
fluids. A set of eigenvalues and eigenfunctions of the generator for this Liouville dynamics is 
identified in a special "stationary representation". This provides a class of exact solutions to the 
Liouville equation that are closely related to hydrodynamics for granular fluids. 
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I. INTRODUCTION 



The properties of simple atomic fluids are typically well-described by classical mechanics 
of point particles with phenomenological pairwise additive potentials. These potentials are 
chosen to capture the quantitative features of the actual quantum mechanical description of 
electronic interactions between pairs of nuclei: strongly repulsive effects at short distances 
and weak attraction at large distances. Except at low temperatures the short range repulsion 
dominates and the classical potentials behave qualitatively as (a/r)~ n where r is the distance 
between a pair, a is a characteristic force range, and n is a large integer (e.g., n = 12 for 
the Lennard- Jones potential). The primary effect of the phenomenological potentials is 
therefore one of excluded volume, dominating all thermodynamic, structural, and transport 
properties. This suggests a final idealization to consider, the hard sphere limit n — > oo, 
where the potential is unbounded for r < a and zero for r > o. The single parameter of the 
hard sphere fluid, a, then can be chosen to give a reasonable quantitative description of real 
fluids |l|. 

The equilibrium statistical mechanics of hard spheres reduces to a determination of con- 
figurations for non-overlapping spheres, a geometrical problem. All thermodynamic and 
structural properties of the hard sphere fluid can be understood simply as a limiting form 
for strongly repulsive potentials. In contrast, the dynamical properties of hard spheres pro- 
vide new problems due to the singular forces and vanishing time scale for pair collisions. 
Consequently, the description of trajectories in the N particle phase space cannot be given 
directly by Newton's or Hamilton's equations Nevertheless, it is clear that 

such trajectories exist and are uniquely defined for the usual initial data. All particles move 
freely until a given pair is at contact, then the momenta of that pair is changed instanta- 
neously according to some chosen collision rule (different for elastic and inelastic collisions, 
as described in the next section). Subsequently the particles continue to move freely until 
the next pair at contact occurs, the collision rule is applied instantaneously to that pair, 
and free streaming is again resumed. The resulting trajectories are unique, deterministic, 
and reversible up to the initial data. Some of the complexities of the crossover from strongly 
repulsive continuous potentials to the hard sphere limit have been described recently jij. 

In the next section a representation of the generators for hard sphere dynamics is recalled 
Such hard sphere generators have been applied with great success over the 
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past thirty years to explore many of the subtleties of transport phenomena in normal fluids 
(e.g., nonanalytic density dependencies, algebraic time decay) (4, (|. More recently, interest 
in hard sphere fluids has been revived as models for granular flows Real granular 

media consist of mesoscopic sized grains (e.g., sand, seeds, pharmaceutical pills) which, when 
activated, behave like complex fluids, i.e., fluids whose constituent particles have non-trivial 
internal degrees of freedom. It turns out that many of the qualitative and quantitative 
differences between normal and granular fluids are captured by a model the hard sphere 
fluid whose constituents undergo inelastic or dissipative collisions. For an overview of recent 
developments see the edited volumes of reference and the new text by Poschel 

Following the definitions of generators, attention is focused in section 3 on their appli- 
cation to the statistical mechanics of granular media. The Liouville equation is described 
and shown to have no stationary solutions (equilibrium) for an isolated system. Neverthe- 
less, a change of variables to accommodate the collisional "cooling" due to the inelastic 
collisions leads to a representation supporting a stationary solution. Next, it is shown that 
the existence of the stationary solution implies certain solutions to the eigenvalue problem 
for the Liouville operator in this representation. In this way a class of special solutions 
to the Liouville equation is identified. It is noted that the eigenvalues are the same as 
those for the exact macroscopic balance equations for mass, energy, and momentum in the 
long wavelength limit. This suggests that the microscopic excitations are the precursors of 
macroscopic hydrodynamics in a granular fluid. 



II. GENERATORS FOR HARD SPHERE DYNAMICS 

The system of interest is a one component fluid of N identical smooth hard disks or 
spheres (mass m, diameter a). The position and velocity coordinates of the fluid particles 
will be denoted by {q^Vj}. The state of the system at time t is completely characterized 
by the positions and velocities of all particles at that time and is represented by a point 
T t = {qi(t), . . . , qjv(t), Vx(t), . . . , Vjv(^)} in the associated 2dN dimensional phase space, 
where d = 2 for hard disks and d = 3 for hard spheres. The dynamics consists of free 
streaming (straight line motion along the direction of the velocity at time t, until any pair 
of particles, say i, j, is in contact. At the contact time the relative velocity gij = Vj — of 
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that pair changes instantaneously according to the collision rule 



Sij = Bij - (1 + a) (5- ■ g {j ) a. (1) 

Here a is a unit vector directed from the center of particle j to the center of particle i 
through the point of contact. The parameter a (the coefficient of normal restitution) is 
chosen a priori in the range < a < 1 and remains fixed for a given system. The value 
a = 1 corresponds to elastic, energy conserving collisions, while a < 1 describes an inelastic 
collision with a corresponding energy loss for the pair. However, the center of mass velocity 
is unchanged so that the total mass and momentum of the pairs are conserved for all values 
of a. Subsequent to the change in relative velocity for the pair i,j the free streaming of 
all particles continues until another pair is at contact, and the corresponding instantaneous 
change in their relative velocities is performed. The sequence of free streaming and binary 
collisions determines uniquely the positions and velocities of the hard particles at time t for 
given initial conditions. 

A. Statistical Mechanics 

The statistical mechanics for a fluid of inelastic hard spheres has been described else- 
where HQ. It is co mpris ed of the d_ just described, a state specified in tenns of a 
probability density p(T), and a set of observables denoted by A(T). The expectation value 
for an observable at time t > for a state p(T) given at t = is defined by 

(A(t);0)= J dTp(T)A(T t ), T t = { qi (*), ■ ■ ■ , q N (t), Vl (i), . . . , v N (t)} (2) 

where A(T,t) = A(T t ), T = T t=0 , and T t = {qi(t), . . . , qjv(0> Vi(£), . . . , v N (t)}. The dy- 
namics can be represented in terms of the generator L defined by 

(A(t);0} = J dT p(T)e tL A(T). (3) 

For continuous potentials the generator is easily recognized from Hamilton's equation as a 
Poisson bracket operation with the corresponding Hamiltonian. However, its identification 
for the discontinuous hard sphere potential is less direct. 

There are two components to the generator, corresponding to the two steps of free stream- 
ing and velocity changes at contact. The first part is the same as for continuous potentials 
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while the second part replaces the contribution from the singular force by a "binary collision 
operator" T(i,j) for each pair i, j 

N -1 N N 

i = Evv, + -EEnM'). (4) 

i=l i=l j^£i 

The binary collision operator T(i,j) for continuous potentials is identified directly from the 
Poisson bracket of Hamilton's equations 

T(i,j) - % = -m- 1 (V w V(|q, - q,|)) ■ (V„ - V WJ ) . (5) 

For hard spheres, the position variables are still continuous functions of time but the mo- 
menta are piecewise constant (in the absence of external forces) and discontinuous. Thus, a 
general phase function has the form 

A(T t ) = 6 (ti — t)A ({qi(t), ■ ■ ■ , qjv(t), vi(0), . . . , Vat(O)}) 

+ E (* " *p) (Vfi - *) ^ ({qi(*)> ■ • • , <lN(t), vi(p), . . • , vjv(p)}) • (6) 
P =i 

Here {t p } are the times for the colliding pairs and {vj(p)} are the velocities in the time 
interval between collisions p and p + 1. The Heaviside theta functions identify the time 
intervals between collisions. Direct differentiation of this form leads to the identification of 
the binary collision operator for hard spheres 

T(i,j) = e(-gij- ■ qoOI&j ■ - CT )(% - !)> ( 7 ) 

where q^ is the relative position vector of the two particles, is the Heaviside step function, 
and bij is a substitution operator 

M(S«) = 4(&<iSy) = A(gy), (8) 

which changes the relative velocity gjj into its scattered value gy, given by Eq. (0). The 
theta function and delta function in (J2J) assure that a collision takes place, i.e. the pair is 
at contact and directed toward each other. Further details of this approach can be found 
in Appendix A of reference [l^ by Lutsko . Alternatively, the same result can be obtained 
through a limiting procedure starting with a continuous potential whose slope increases 
without bound. 
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An alternative equivalent representation of the dynamics is to transfer it from the ob- 
servable A(T) to the state p(T) by the definition 

/ dT p(T)e tL A(T) =jdT (e- tT p(T)) A(T). (9) 

The representation in terms of a dynamical state is referred to as Liouville dynamics. Implicit 
in the analysis of the previous paragraph for hard spheres is the restriction of the phase 
space to non-overlapping configurations. This is assured when the generator L is used in the 
context of averages such as (jSJ) since all acceptable probability densities p(T) must include a 
singular factor excluding the domain of any overlapping pair. However, the right side of (jHJ) 
no longer has that restriction and consequently the generator for Liouville dynamics is not 
the same as that for observables (as in the case of continuous potentials). Instead, direct 
analysis of (0) leads to the result 

N 1 N N 

z = £vi-v EE^'). (io) 

i=l i=l j^i 

with the new binary collision operator 

T-(i,j) = <f(g y - <7)|gjj • qtf|(6(gu • q^oT 2 ^ 1 - 6(-g ii • %■))■ (H) 

Here 6" 1 is the inverse of the operator 6^- in (jHJ) 

1 + a 

KjSij = Sij (o- ■ Sij) or- (12) 

a 

In summary, the problems presented by the singular forces for a fluid of hard spheres are 
resolved if Hamilton's equations for observables are replaced by 

(d t -L)A(T,t) = 0, (13) 

and the Liouville equation for probability densities is replaced by 

(d t + L)p(T,t)=0, (14) 

with the respective generators given by (J3J and (|10p. An important observation empha- 
sized by Lutsko [12| is that the form of the generator L and corresponding binary collision 
operator T(i,j) does not depend on the details of the collision rule defining the operator 
bij. In particular, the result applies for both elastic and inelastic collisions. In contrast, 
the generator for Liouville dynamics is obtained by a change of variables that introduces 
the Jacobian of the transformation between the variables gjj and ftygy. Hence it depends 
explicitly on the collision rule, or for granular media on the restitution coefficient a. 
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III. LIOUVILLE EQUATION FOR GRANULAR FLUIDS 



In the remainder of this presentation attention will be focused on the granular fluids case 
a < 1, for an isolated system. In contrast to normal fluids, there is no stationary solution to 
the Liouville equation for an isolated system. This follows by calculating from it the average 
"thermal" speed of the particles v (t), defined in terms of the average kinetic energy 

■fi4^1™»| (15) 

where m is the mass, and the second equality defines the average . Using either representa- 
tion or ([T3J) the average speed is found to be monotonically decreasing 

d t \nv(t) = ~((t), (16) 

where ((t) > is the "cooling" rate due to inelastic collisions, 

AT r 

C (t) = (1 - « 2 )^y J rfr(g 12 -q 12 ) 3 6(g 12 • q 12 )%i2 - <T)p(T, t). (17) 

This shows that there is no "approach to equilibrium" for a granular fluid since there is no 
such stationary equilibrium state. Nevertheless, it is expected (on both theoretical grounds 
and computer simulation results) that there is a universal state that is approached for a 
wide range of initial preparations. 

This universal state, the homogeneous cooling solution (HCS), is special in the sense 
that it is spatially homogeneous (translationally invariant) and all of its time dependence 
occurs through the average speed v(t). If there are no other externally imposed energy scales, 
dimensional analysis requires that the dependence on v(t) occurs only through normalization 
and through a scaling of the velocities, 

*- = (**-(«))^A.({f.^})- («) 

Here Vh cs {t) is the thermal speed evaluated for the HCS, u is the constant average velocity 
of the system and £ is a constant characteristic length (e.g., the mean free path). In this 
state, collisional cooling due to inelastic collisions amounts only to a monotonic decrease 
in all velocities. Use of ()18j) in (JT7j) shows that the cooling rate in the HCS has the time 
dependence 

a cs (t) = ^#c cs , (19) 



where Q cs is a dimensionless constant. Then (frfijl can be solved directly to get the explicit 
time dependence of Vh cs {t) 

v hcs (t) = v hcs (0) (l + ^ffW) 1 - (20) 

Also given is the limiting behavior at long times, showing that the scaling becomes indepen- 
dent of the initial conditions. 



A. Stationary Representation 

More generally, collisional cooling affects both the form of the distribution function as 
well as scaling the velocities. It is useful to account for the latter by making a time dependent 
change of variables to obtain a dimensionless form of the Liouville equation 

q,=y, ^=^W ds = — dt ' ( 21 ) 

P*({q*,va, s ) = (^(t)) 37V p(r,t). (22) 

Here u(t) is a characteristic velocity whose form is to be chosen below. In terms of these 
variables the dimensionless Liouville equation becomes 

uj (t) N — 
d sP * - l^rfr E Vv* • (V*p*) + Vp* = 0, (23) 
u [t) i=l 



with the time independent dimensionless generator 



(24) 



Next, uj(t) is chosen to make the coefficients of the Liouville equation independent of s 

where uj* is an arbitrary dimensionless constant taken here to be unity. The solution to this 
equation is 

u{t)=u,(0)(l + ^tj , (26) 
and consequently the new time scale is 

+ (27) 



The final form for the dimensionless Liouville equation is now 

d,p*+Tp* = 0, (28) 

with the modified generator 

N 

Tp* = Tp* + J2 v v * • (V t V) (29) 



1=1 



This time dependent change of variables was first suggested by Lutsko [ljj , motivated by 
the objective of molecular dynamics simulation of the HCS. Direct simulation in the original 
variables is difficult at long times since the velocities monotonically scale toward zero. In 
contrast, the corresponding solution to ()28|) is now a stationary solution p^ given by 

Tp* = 0. (30) 

Consequently, the corresponding thermal speed for the new variables approaches a constant 

/£l * f «. n \ _ 2 - 



r *2 



2 r N 
-> — dr>* (T) £ vf = ^ 2 (oo) (31) 

Since both Vh cs (t) and u;(i) decay as t _1 at long times it follows that the dimensionless 
cooling rate Q cs is simply related to v*(oc) 

2 

Chcs = v^(ooj' (32) 

Since (|2*Hj) supports a stationary state it is referred to here as the stationary representation 
of the Liouville equation. 

B. Special Solutions to the Liouville Equation 

The existence of the stationary state p^ of the generator C* given by (j3*Uj) has some in- 
teresting consequences. The stationary solution is assumed known. In the original variables 
it satisfies the equations 

p (T,t) = (iu (t))- 3N pl ({a*, ^^}) (33) 
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LPo + ^rY, V Vi ■ ((v, - m) po) = 0. (34) 
Differentiating (jHjj) with respect to u; (£) at constant {q*,Vj} gives 
— cjM N N 

Ldupo + -y- J2 V Vl ■ ((Vi - Ui) d wPQ ) = -f l Y, V v. ■ ((v< - u) Po). (35) 

*- i=l i=l 

But it also follows from (J33|) that 



N 

d^Po = -lj^Y V Vl ■ ((vi - Uj)p )- (36) 

i=l 

Substitution of into (J35J) and transformation back to dimensionless variables gives the 
result 

N 

£*^(i) = ^(i) = £ Vv . . (y* p *). (37) 

i=l 

Thus, the stationary state implies directly one solution to the eigenvalue problem for C . 

Repeating this analysis, but differentiating with respect to the total number of particles 
iV and the flow velocity u gives the additional eigenvalues and eigenvectors 

ZV 2 ) = 0, ¥ 2 > = d N p* (38) 

N 



1=1 



The eigenvalue —1 is three fold degenerate since the equation holds for each component of 
the vector Hf( 3 \ 

These eigenvalues and eigenvectors allow construction of a special class of exact solutions 
to the Liouville equation (|28j) 

p* (r, s) = e - Ts [pi (r*) + ci^ (1) (r*) + c 2 ¥ 2) (r*) + c 3 • * (3) (r*)) 

= p* (r*) + ci^ (1) (r*) e~ s + c 2 ^ (2) (r*) + c 3 • * (3) (r*) e s (40) 

where {q} are arbitrary constants. Note that normalization is preserved since the integral 
over phase space of all eigenfunctions vanish. In terms of the original variables this class of 
solutions can also be written 

p (r, s) = po (T) + £ e~ x " s c n d ynPo (r) (41) 

n=l 

with 

A n — (1,0, -1, -1, -1) , y n — (u, N, u) . (42) 
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C. Relationship to Hydrodynamics 



To provide some interpretation to the dynamical effects represented by these eigenfunc- 
tions and eigenvalues it is useful to consider the exact macroscopic balance equations for 
density n, momentum density mnu, and energy density ^nksT. These follow directly by 
using either (|T3j) or (|TU|) to calculate the averages of the corresponding local microscopic 
densities, with the results 

D t n + nV • u = 0, (43) 

D t Ui + {mn)~ 1 d i p + (mn)~ 1 djPij = 0, (44) 

(A + C) T + J-pV • u + £ (PijdjUi + V ■ q) = 0, (45) 

where D t = dt + u • V is the material derivative, T(r, t) is the temperature, u(r, t) is the 
flow velocity, q is the heat flux, and is irreversible part of the momentum flux. Linearize 
these equations about the HCS solution {n hcs = constant, T hcs {t) = mvf lcs {t)/2, u hcs = 0), 
and introduce the dimensionless variables 

x * 6n vr* 6T x * 6n (Aa\ 
Sn = , ST = 5u = — . 46 

n hcs T hcs (t) v hcs (t) 

Then in the long wavelength limit, neglecting all spatial gradients, these equations become 

id, + 0) 5n* w 0. (47) 

(d s + 1)( 5T* + d ^ n(hcs 5n*) « 0. (48) 
V d\nn hcs J 

(d s -l)5U*^0. (49) 

Therefore, it is seen that excitations obtained for the Liouville equation are the same dy- 
namics as the long wavelength limit for the macroscopic balance equations. 

More generally, it is expected that spatially inhomogeneous states for a granular 
gas can be described by hydrodynamics. Such hydrodynamic equations follow from the 
above macroscopic balance equations with constitutive relations expressing the heat flux 
and momentum flux in terms of spatial gradients of n, u, and T. The resulting linearized 
equations define the hydrodynamic modes for a granular system, which necessarily behave 
as (J47)) - (J48j) at long wavelengths. This shows that the eigenfunctions and eigenvalues of 
the Liouville equation found here can be interpreted as the long wavelength microscopic 
precursors of macroscopic hydrodynamics. This point is illustrated in more detail elsewhere 
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IV. CONCLUSION 



The objectives of this brief presentation have been two fold. The first was an introduction 
to the means to describe the dynamics of a hard sphere fluid, in spite of the singular forces 
and failure of standard Newtonian and Hamiltonian formalisms. In fact, this is a special 
case of a class of problems involving piecewise continuous or discontinuous potentials (e.g., 
the finite step potential, hard spheres with attractive square well)j3|. In each case the effect 
of the singular force can be replaced by an associated binary collision operator to define a 
generator for the dynamics. 

The second objective was to report some new results for the Liouville dynamics of inelas- 
tic hard spheres, an idealized model for granular fluids. The inelasticity implies a monotonic 
loss of energy for an isolated system. However, by a change of variables to accommodate the 
average decrease in each particle's speed the Liouville equation is given a representation that 
supports a stationary state. The existence of this stationary state implies certain properties 
of the associated generator for the dynamics. For systems with elastic collisions this would 
imply certain invariants. For inelastic collisions the eigenfunctions identified here are no 
longer invariants, but instead have a simple dynamics directly related to the collisional cool- 
ing rate. This can be seen by choosing as the scaling function u>(t) — > Vhcs{t)- In this case the 
eigenvalues become A n < — > [hChcsj 0; — ^Chcsi ~ \Chcai ~ |Cftcs) • It was then noted that 
the eigenvalues of the Liouville equation are the same as those of the macroscopic balance 
equations for average number density, energy density, and flow velocity in the long wave- 
length limit. Hence they are long wavelength hydrodynamic modes, whenever a closed set 
of hydrodynamic equations apply. This has been exploited recently to define hydrodynamic 
response functions for a granular fluid and to derive corresponding exact Green-Kubo expres- 
sions for the transport coefficients 14] . This identification of eigenvalues and eigenvectors of 
the Liouville operator is a generalization of a similar analysis at the level of the Boltzmann 
kinetic theory The connection to hydrodynamics has been established clearly in that 
context. 
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